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ABSTRACT 

This  report  examines  the  problem  of  estimating  the  location*  velocity  and 
manoeuvre  of  a  manoeuvring  target  (from  range,  bearing  and  pose  informa¬ 
tion).  The  report  considers  a  non-linear  modeling  technique  in  which  the  target 
is  represented  as  a  hybrid  system  (a  combination  of  discrete  and  continuous  val¬ 
ued  states)  and  considers  new  associated  approaches  such  as  the  polymorphic 
estimator. 

Although  simulation  studies  were  performed,  the  polymorphic  estimator 
had  serious  numeric  problems  that  suggested  the  estimator  should  not  be  used 
until  the  approach  is  refined.  This  report  is  intended  to  facilitate  further 
discussion  and  development  of  the  approach  (if  deemed  necessary);  hence,  it 
includes  details  of  the  assumptions  made  and  the  Matlab™  implementation. 
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Estimation  of  Manoeuvring  Targets  using  Hybrid  Filters 

EXECUTIVE  SUMMARY 

Modern  guided  weapons  are  frequently  required  to  operate  in  a  complex  environment  that 
often  involve  highly  complicated  behaviours.  In  these  types  of  engagements,  assumptions 
of  linearity  no  longer  hold  and  model  uncertainties  in  the  form  of  unmeasured  aerody¬ 
namic  coefficients  and  complex  non-linear  aerodynamics  are  common.  The  target  filter 
is  an  important  sub-system  of  any  guidance  loop  that  estimates  the  required  target  and 
engagement  information.  To  improve  the  performance  of  this  target  filter  in  challenging 
environments  involving  manoeuvring  target,  a  full  understanding  of  any  non-linearities 
present  is  required. 

The  aim  of  this  report  is  to  investigate  the  manoeuvring  target  filtering  problem*  to  exam¬ 
ine  the  importance  of  mode  measurements  and  to  examine  a  particular  filtering  approach. 
A  review  of  existing  filtering  results  is  provided  before  introducing  a  non-linear  filtering 
approach  known  as  hybrid  filtering.  Three  possible  filtering  approaches  are  examined  in 
simulation  studies:  the  extended  Kalman  filter,  the  interacting  multiple  model  filter  and 
a  new  hybrid  filtering  approach.  The  simulation  studies  suggest  that  mode  measurements 
may  improve  target  filtering  performance  but  the  studies  do  not  support  the  use  of  the 
examined  hybrid  filtering  approach.  Some  refinement  of  this  hybrid  filtering  approach  is 
required. 

An  improved  understanding  of  filtering  techniques  is  required  to  aid  support  of  present 
upgrade  programs  involving  the  guidance  loops  of  new  air-to-air  and  standoff  missile  sys¬ 
tems.  This  understanding  is  necessary  for  the  support  of  future  weapon  procurement  and 
upgrade  programs. 
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1  Introduction 

Precision  guidance  of  weapon  systems  is  a  computationally  and  conceptually  demanding 
problem  [2].  Historically,  due  to  real-time  computing  constraints,  major  approximations 
in  the  control  design  process  have  been  necessary.  Recent  advances  in  missile  sub-systems 
mean  that  modern  guided  weapons  have  significantly  improved  computational  capacity 
and  hence  various  modelling  approximations  are  being  reconsidered  [2]. 

One  sub-problem  of  the  missile  guidance  problem  is  estimation  of  the  target  position  and 
velocity  (and  other  quantities)  from  measurements  such  as  range  and  bearing  to  the  target. 
A  common  approach  involves  developing  a  system  model  of  the  relationship  between  the 
target  state  (position  and  velocity  say)  and  the  measurements  available.  Once  a  system 
model  has  been  obtained,  the  target  estimation  problem  can  then  be  posed  as  a  model 
based  filtering  (or  optimal  filtering)  problem. 

Across  many  fields  of  study,  one  of  the  more  famous  and  commonly  used  system  models 
describes  the  relationship  between  the  states  of  a  system  and  system  measurements  as  a 
linear  Gauss-Markov  system  with  Gaussian  noises.  This  model  assumes  linear  dynamic  be¬ 
haviour  of  the  internal  system  state  (with  perhaps  an  additive  Gaussian  noise  disturbance) 
and  measurements  that  are  noisy  linear  functions  of  the  state. 

This  linear  Gauss-Markov  system  assumption  is  popular  because,  although  most  systems 
have  no  finite  dimensional  optimal  filter,  the  Kalman  filter  has  been  shown  to  be  the 
optimal  filter  for  such  a  system  [1,  5].  In  this  context,  optimality  is  in  the  minimum 
mean  squares  sense  and  a  filter  is  finite  dimensional  If  it  can  be  implemented  using  a  finite 
number  of  statistics  which  can  be  calculated  using  a  finite  number  of  recursions.  Because 
the  Kalman  filter  is  a  finite  dimensional  optimal  filter  it  has  been  applied  to  a  large  variety 
of  filtering  problems. 

In  a  typical  Interceptor-target  engagement,  measurements  are  the  relative  range  and  bear¬ 
ing  to  the  target  (and  perhaps  pose  information).  These  measurements  are  non-linearly 
related  to  the  relative  position  and  velocity  of  the  target  in  cartesian  co-ordinates.  Because 
of  the  non-linearity  in  the  measurement  process,  the  Kalman  filter  Is  not  appropriate  for 
this  problem  (even  under  the  assumption  that  the  state  dynamics  of  the  target  are  Gauss- 
linear).  Further,  if  the  target  is  performing  manoeuvres  then  a  Gauss-linear  model  of  the 
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target  dynamics  is  no  longer  appropriate. 

For  general  non-linear  problems,  when  a  finite-dimensional  optimal  filter  is  not  possible, 
sub-optimal  numeric  or  approximate  approaches  must  be  used.  The  simplest  approach 
is  to  use  an  extension  of  the  Kalman  filter  known  as  the  extended  Kalman  filter  (EKF) 
[5,  1].  The  EKF  involves  linearisation  of  the  non-linear  model  about  the  current  operation 
point.  The  EKF  approach,  although  appealing  due  to  its  similarity  to  the  Kalman  filter, 
does  not  perform  well  in  many  non-linear  filtering  problems. 

Some  non-linear  systems  can  be  well  represented  by  a  non-linear  model  known  as  the 
hybrid  system  model  (a  model  that  contains  a  mixture  of  continuous  and  discrete  valued 
states).  Close  to  optimal  finite  dimensional  filters  have  been  designed  that  perform  better 
than  EKF  approaches  on  these  systems. 

In  this  report  we  compare  three  non-linear  filtering  approaches  to  the  problem  of  estimat¬ 
ing  the  state  of  a  manoeuvring  target.  The  three  approaches  are  a  simple  EKF  approach 
(included  as  a  benchmark),  the  interacting  multiple  model  (IMM)  filter  [9],  and  the  poly¬ 
morphic  estimate  (PME)  [18].  The  IMM  and  PME  are  filtering  approaches  designed  on 
an  approximate  hybrid  model  representation  of  the  manoeuvring  target.  An  important 
issue  is  whether  these  hybrid  system  models  are  realistic  representations  of  manoeuvring 
target  dynamics. 

The  key  aim  of  this  report  is  to  examine  the  performance  of  these  three  filtering  approaches 
in  a  manoeuvring  target  filtering  problem  to  evaluate  the  apparent  advantages  of  a  hybrid 
system  filtering  approach.  The  second  aim  of  this  report  is  to  examine  these  filters  as  part 
of  a  missile  guidance  loop.  A  third  aim  is  to  examine  the  influence  of  pose  information  on 
a  missile  guidance  loop. 

The  report  is  structured  as  follows:  In  Section  2,  we  introduce  several  possible  dynamic 
models  and  then  introduce  the  target  filtering  problem.  In  Section  3,  the  extended  Kalman 
filter  solution  to  the  problems  is  presented.  In  Section  4,  the  interacting  multiple  model 
filter  and  the  polymorphic  estimator  are  presented.  In  Section  5,  some  implementation 
issues  are  discussed.  In  Section  6,  simulation  results  are  presented.  Finally,  in  Section  7, 
some  conclusions  are  presented. 
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2  Target  Filtering  Problem 

In  this  section  we  present  the  manoeuvring  target  tracking  problem.  Many  formulations  of 
the  target  tracking  problem  are  commonly  used.  The  formulation  presented  below  allows 
several  different  types  of  dynamic  models  to  be  considered  in  a  similar  formulation.  We 
first  introduce  the  filtering  problem,  then  introduce  two  target  models. 

2.1  The  Filtering  Problem 

The  manoeuvring  target  filtering  problem  stated  in  the  broadest  terms  is  to  determine  in¬ 
formation  about  the  state,  usually  denoted  xk^  and  perhaps  the  target  manoeuvre,  denoted 
here  from  measurements  up  until  time  k.  In  the  context  of  the  target  tracking  prob¬ 
lem,  we  are  usually  interested  in  estimating  the  position,  velocity  and  perhaps  manoeuvre 
of  the  target  from  range  and  bearing  measurements. 

For  the  purpose  of  this  report  we  are  going  to  consider  filtering  and  estimation  from  a 
model  based  perspective  where  estimation  is  according  to  a  conditional  mean  criteria  (we 
will  limit  our  interest  to  statistics  of  the  state  such  as  the  mean  and  variance). 

An  exact  study  of  the  properties  of  the  presented  filters  is  beyond  the  scope  of  this  re¬ 
port.  Various  practical  approaches  are  considered  in  the  following  sections  but  none  of 
these  filters  are  optimal.  We  will  examine  the  proposed  approaches  according  to  their 
performance  in  a  target  tracking  problem. 

2.2  Continuous  State  Model 

Consider  the  following  non-linear  model  for  engagement  with  a  manoeuvring  target: 

=  ak(xk,  u[ ,  uf,  vk) 

% ik  =  ck(xk,wk)  (2.1) 

where  xk  £  RN  is  a  state  describing  both  the  target  and  interceptor  dynamics  (or  perhaps 
the  relative  dynamics  between  the  two),  £  U{xk)  is  the  manoeuvre  performed  by  the 
target  which  is  an  unmeasured  input,  uk  is  the  manoeuvre  performed  by  the  interceptor 
which  is  assumed  to  be  known,  yk  £  RM  is  the  observation,  vk  is  a  process  noise  process, 
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and  Wk  is  a  measurement  noise  process.  The  set  U(xk)  is  the  complete  set  of  target 
manoeuvres  that  can  be  performed  from  state  xk. 

This  model  is  generic  enough  to  include  general  target  manoeuvres,  general  measurement 
processes,  and  complicated  aerodynamics.  The  problem  considered  in  this  report  is  esti¬ 
mation  of  the  state,  x *,  from  the  measurements  y k.  The  general  filtering  problem  for  the 
above  non-linear  system  with  unmeasured  uj  has  not  been  solved  in  analytic  form  (there 
are  several  numeric  approaches  such  as  the  particle  filter  which  are  outside  the  scope  of 
this  report  and  not  considered  here). 

In  some  situations  it  is  reasonable  to  use  a  more  restrictive  stochastic  process  to  describe 
v%  and  to  assume  linear  state  dynamics  as  follows: 

xk+i  =  Axk  +  BTuTk  +  BJurk  +  vk 
Vi k  =  c*(a*)  +  n*  (2.2) 

where  A  G  RNxN ,  vk  and  nk  are  Gaussian  noise  processes,  and  i%  is  a  stochastic  process 
dependent  on  xk  and  v%_v 

This  filtering  problem  is  also  difficult  and  no  analytic  solution  exists.  The  model  (2.2)  is 
appropriate  for  an  extended  Kalman  filter  solution  to  the  target  estimation  problem.  We 
will  discuss  this  approach  later  but  first  we  introduce  further  restrictions  on  the  input  that 
lead  to  a  hybrid  system  model. 


2.3  Hybrid  System  Model 

In  this  subsection  we  consider  a  further  restriction  on  uj  that  leads  to  a  hybrid  system 
model  and  leads  to  two  alternative  filtering  approaches.  Consider  the  case  that  U(xk)  is 
restricted  to  be  a  finite  set  (whose  elements  correspond  to  distinct  possible  target  acceler¬ 
ations)  and  further  that  v%  is  a  Markov  process.  With  this  finite  set  restriction,  assuming 
linear  state  dynamics,  the  following  model  is  obtained: 

1 4  is  a  first  order  Markov  chain 

xk+i  =  A(uJ)xk  +  Brurk  +  vk 

Vk  =  ck(xk)  +  wk.  (2.3) 
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where  is  a  finite  proems  that  takes  values  from  a  discrete  set  (of  size  Nv)  whose  values 
approximate  tlie  behaviour  of  the  target.  Here  each  A(.)  indexed  by  uj^  describes  the 
target  dynamics  for  each  of  the  possible  target  manoeuvres.  For  example,  A(0)  —  A  where 
0  indicates  the  target  is  not  performing  a  manoeuvre  etc.  It  should  be  noted  that  this 
target  acceleration  model,  (2.3),  only  holds  when  the  possible  target  control  actions  are 
restricted  to  perpendicular  accelerations. 

This  system  is  considered  a  hybrid  system  model  because  it  contains  a  mixture  of  contin¬ 
uous  valued  and  discrete  valued  states  (xk  and  ujT  respectively). 

The  motivation  in  considering  a  hybrid  system  approximation  for  the  continuous  system 
(2.1)  stems  from  the  knowledge  that  “close  to  optimal55  filtering  solutions  for  hybrid  sys¬ 
tems  exist.  The  meaning  of  “close  to  optimal55  will  be  discussed  in  a  later  section.  A 
key  question  is  whether  a  “close  to  optimal55  solution  on  the  hybrid  system  approximation 
achieves  better  performance  than  an  approximate  solution  based  on  a  continuous  system 
model.  This  issue  is  investigated  in  Section  6  via  simulation  studies. 

The  system  model  (2.3)  is  appropriate  for  both  the  interacting  multiple  model  and  the 
polymorphic  estimator  (both  are  discussed  in  Section  4). 

3  The  Extended  Kalman  Filter 

The  Kalman  filtering  is  the  optimal  filter  for  a  linear  Gauss-Markov  system.  System  (2.2) 
is  linear  Gauss-Markov  when  A(uk)  —  A/.,  Ck{xj-)  =  CkXk  and  Vk,Wk  are  independent 
sequences  of  Gaussian  noise.  The  Kalman  filter  is  optimal  in  the  sense  that  it  produces 
estimates,  that  minimise 

E[(xk  -  £fc|fc)(:r*  -  *  *  • ,  Vk\- 

The  extended  Kalman  filter  extends  the  concept  of  the  Kalman  filter  to  a  non-linear 
system  model  via  linearisation.  The  extended  Kalman  filter  (analogous  to  the  Kalman 
filter)  calculates  a  state  estimate,  a  covariance  matrix,  a  priori  state  estimate  and  a  priori 
covariance  matrix  at  each  time  instant  (£*|jt,P*|fc,  xk\k-\  and  Pk\k-i  respectively). 

To  develop  the  extended  Kalman  filter  for  this  problem  we  first  assume  that  u[  is  available 
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£*|k-i  =  0*-l(^t-l|fc-l,  tfc  i  uib  vk) 

■Pfc|fc-i  =  A.k-\Pk-\\k-\A'k_x  +  Q*k 

Kk  =  Pk\k-iC'k[ckPklk_lC'k  +  R]\~l 

=  1  "I"  K k  ^J//s  ck(&k\k— l)] 

Pk\k  =  Pk\k-i  -  KkCkPk\k-\  (3.3) 


Given  a  sequence  of  measurements,  t/*,  the  extended  Kalman  filter  provides  a  sequence  of 
estimates,  and  an  estimate  of  error  covariance  Although  in  non-linear  estimation 
problems  higher  order  moments  can  be  significant,  the  EKF  only  keeps  track  of  1st  and 
2nd  moment  information  {xk\k-uPk\k-i)- 

The  above  implementation  assumes  that  an  estimate  of  the  target  manoeuvre  is  available. 
Often  a  target  manoeuvre  estimate  will  not  be  available  and  a  common  strategy  is  to 
then  assume  the  target  is  not  manoeuvring.  This  sort  of  approximation  about  the  target 
manoeuvre  significantly  affects  the  performance  of  the  filter  and  is  the  prime  motivation 
for  considering  the  two  filters  described  in  the  next  section. 

Another  possible  alternative  is  to  include  the  target  manoeuvre  as  part  of  the  system  state. 
However,  when  formulated  this  way  the  system  model  becomes  highly  non-linear  and  the 
extended  Kalman  filter  itself  (which  is  based  on  linearisations)  is  likely  to  be  unstable 
when  initialised  poorly. 
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4  Hybrid  System  Filtering 

In  this  section  two  separate  approaches  to  the  filtering  problem  for  hybrid  systems  are 
discussed.  The  first  filtering  approach  considered  here  is  the  interacting  multiple  model 
filter  which  is  based  on  the  concept  of  a  running  bank  of  filters  (one  for  each  of  the  possible 
manoeuvre  modes)  and  then  combine  the  output  of  each  filter  to  obtain  an  estimate  of  the 
target  state.  The  second  filter  considered  is  the  Polymorphic  Estimators  which  is  based 
on  a  sub-optimal  approximation  to  the  optimal  solution  of  the  hybrid  system  filtering 
problem. 


4 .1  Interacting  Multiple  Model  Filter 

The  interacting  multiple  model  filter  is  a  filtering  approach  (based  on  the  model  (2.3)) 
that  involves  running  a  sub-filter  for  each  of  the  distinct  manoeuvres  modelled  [9].  At 
each  time  step,  each  individual  sub-filter  uses  any  new  measurements,  mode  probability 
information  and  previous  sub-filter  outputs  to  generate  a  new  estimate  based  on  the  sub¬ 
filter  ?s  assumed  manoeuvre.  Then  the  output  all  these  sub-filters  is  combined  according 
to  certain  probability  information  to  produce  a  new  current  overall  target  estimate.  This 
is  explained  in  more  detail  below. 

Let  and  P|_1|^_1  be  the  outputs  of  the  ith  sub-filter  at  time  k  —  1  and  let  A& 

be  the  transition  probability  from  manoeuvre  i  to  manoeuvre  j.  Also,  let  denote  the 
estimated  probability  of  being  in  mode  j  at  time  k  —  1. 

Then 


a 

4V. 


pjo 

■Nfc-llfc-l 


N+U 

i-l 
Nv 


where  and  P^_  are  the  inputs  to  the  jih.  sub-filter  at  the  Arth  cycle  used 

below,  is  the  probability  that  a  particular  transition  of  manoeuvre  occured,  and  cP 
is  a  normalisation  constant  (from  p^_i  over  i). 
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Then  the  Nu  sub-filters  iterate  as  follows 

4\k-i  =  o*-i(4°-i|fc-i>ufc  =j.  «*.«*) 

H\k-i  = 

K  - 

”  4|fc-l  [yk  “  cA:(^i|fc-l)] 

^|*  =  Pi\k-,-HCkPi^  (4.2) 

where  C*  was  defined  in  the  last  section  and  j4*(j)  =  =  i). 

The  mode  probability  can  be  calculated  as 

H  =  P  {Vkl^k  =  J^k-i\k-VPk-l\k-l) 

Nu 

t*k  =  cAiE^X- 1  (4.3) 

t=l 

where  c  is  a  normalisation  constant. 

The  overall  estimate  of  the  overall  scheme  at  any  time  is  given  by 

Nv 

4*i  k  =  E^i  A 

3=1 

Nu 

Pk\k  =  E^i(pfc|fe  +  (4|fc-**|fc)(4|fc-^|fc)')  (4-4) 

j- 1 

Note  that  the  IMM  filter  keeps  track  of  only  1st  and  2nd  order  moment  quantities 


4.2  Polymorphic  Estimator 

Unlike  the  EKF  and  IMM  filters  presented  above,  the  polymorphic  estimator  keeps  track 
of  several  higher  order  statistics  to  improve  the  quality  of  state  estimates  for  the  hybrid 
system  (2.3).  By  studying  the  optimal  filtering  problem  for  (2.3)  and  making  approxima¬ 
tions  in  equations  involving  4th  and  5th  order  moments,  filters  with  improved  performance 
can  be  obtained.  We  provide  no  other  details  about  the  development  of  the  PME  and  the 
reader  is  directed  to  [18]  for  any  more  information.  Here  we  simply  present  the  filter. 

To  simplify  the  presentation,  we  consider  the  PME  as  a  two  step  process:  a  time  update 
step,  and  a  measurement  update  step. 
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Time  Update 

la  between  measurements,  the  statistics  of  the  hybrid  system  evolve  according  to  the 
underlying  dynamics  of  the  system.  For  this  reason,  the  time  update  step  of  the  PME  is 
typically  done  via  numeric  integration  of  the  following  differential  equation.  In  general, 
the  evolution  of  these  equation  can  not  be  implemented  as  difference  equations. 

P't+st  =  Ay.t 
d  Nu 

jt&t  =  Y,At(i)  (%|+pw) 

i—l 

d  Nu 

dlpx/i  =  Px/iAc  +  23  At(i)  (&tPf$  +  Px^ i  +  Pxfififj 
d  NV  ' 

~dl^>xx  ~  5Z  (f'tPpz  +  PxXxy. *  +  P txfJ't)  +  R{i)n\ 
d  Nv 

-^jPxxii™  —  23  At(i)P(x/J)xiim  +  Plzti')xftrnAt{iy  +  PXXfliQtm  (4-5) 

i—l 

where  is  the  it h  row  of  A  and  A^  is  the  iih  column  of  A.  Note  that 

Pm  -  diagfat)  -  lHl4 

—  PfJLX 

Pxpti*  ~  (ei  —  fit)  + 

—  Pzxflifi™  +  fAtPxXflm  +  StPpixfl™  ~~  Pxp  (Px™^ 

PpL iXflm  ~  PxftpF1  an<^ 

p  .  =  /  +  Pxxfi171  -  fimP tXpi  -  {AP xxfi™  ~  PxxfAfim  if  i  =  m  (  . 

\  -^PXXili  -  ^PXXfim  -  otherwise  (46) 

where  diag{X)  is  the  diagonal  matrix  with  X  on  its  diagonal. 


Measurement  Update 


If  either  range  and  bearing  or  pose  information  becomes  available,  then  the  state  estimate 
of  the  hybrid  system  filter  can  be  updated  as  shown  in  the  following  two  sub-sections. 
Note  that  if  both  types  of  measurements  become  available  at  the  same  time  instant  then 
the  update  steps  can  be  performed  sequentially.  Although  the  following  equations  are 
generally  implemented  as  discrete  time  recursions,  to  ensure  conformity  with  the  time 
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update  equations,  we  have  expressed  the  updates  as  changes  from  xt  to  xf  etc.  at  the 
discrete  time  instants  that  measurements  occur. 

Mode  Measurements  Assuming  that  there  are  measurements  of  the  target  pose  y™, 
this  information  can  be  incorporated  into  the  filter  estimates.  Here  y™  E  {ei,...}  is  a 
discrete  value  representation  of  the  target  pose  (for  example,  the  target  orientation  in  the 
yaw-plane). 


y?  =  Amy? 

pt  =  Pi  * 

An  =  nt  -  nt 

xt  =  x;  +  PxtlyT\ 

Ax  =  P tfiVt 

Nv 

PxpL  —  Px^i  —  Ax  A  y!  +  PX^VT  (0 

Nv 

Pxx  =  Pxx  -  AxAx  +  Y^Pxx^yFi*) 

t=i 


where  the  ij th  element  of  Am  is  the  probability  of  being  in  mode  i  given  that  the  target 
pose  is  j,  iff1  is  a  pseudo  measurement  of  target  mode,  and  *  is  the  operation  of  element- 
by-element  multiplication. 

Range/bearing  Measurements  Range  and  bearing  information,  yt>  can  be  used  to 
update  the  target  state  information  as  follows 

7*  =  PXxC'(CPxxC' +  R)-1 

xt  =  K  +7 xVt 

P%n  —  Pxn  —  'IfxCP xfi 

Pxx  =  Pxx  7x{C  PxxP  "i"  R)j’x  (4-®) 

4.3  Other  Filtering  Approaches 

There  are  many  other  non-linear  filtering  approaches  that  may  be  considered  appropriate 
for  this  problem.  The  robust  Kalman  filtering  approach  may  be  considered  an  alternative 
to  the  extended  Kalman  filtering  approach  [16].  More  success  can  be  obtained  from  non¬ 
model  based  approaches  such  as  Monte  Carlo  approaches  [15]  and  the  particle  or  auxiliary 
particle  filter  approach  [17]. 
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Also  of  interest  is  the  collection  of  approaches  to  the  hybrid  filtering  problem.  There  are 
many  sub-optimal  hybrid  filtering  approaches  available  and  this  paper  only  had  time  to 
examine  one.  Some  of  the  other  hybrid  filtering  approaches  many  be  worth  considering. 

None  of  these  alternative  approaches  were  examined  in  the  context  of  this  simulation  study, 
but  previous  work  suggests  that  the  particle  filter  approach  may  offer  some  advantages 
(once  the  technique  has  become  well  developed).  Further  comment  is  outside  the  scope  of 
this  report. 


5  Implementation  Issues 

This  section  discusses  some  of  the  implementation  issues  associated  with  the  above  filters. 
In  the  first  sub  section  we  will  discuss  some  issues  directly.  Further  information  about  the 
Matlab™  implementation  and  the  syntax  and  header  files  can  be  found  in  Appendix  A. 

5.1  Direct  Issues 

Implementation  of  the  extended  Kalman  filter  has  been  covered  by  many  other  authors  (See 
[5]  for  example),  and  implementation  issues  with  respect  to  the  IMM  are  well  understood 
[9]  so  this  section  will  concentrate  on  issues  related  to  the  Polymorphic  Estimator. 

Numeric  Stability- 

In  our  studies,  the  Polymorphic  Estimator  was  found  to  have  serious  numeric  stability 
problems.  The  most  significant  of  these  problems  was  the  instability  in  the  time  update 
step  of  the  filter.  The  time  update  step  requires  numeric  integration  of  a  non-linear 
system  of  equations.  In  our  implementation  of  the  filter,  a  Euler  approximation  was  used 
for  integration  purposes. 

It  was  found  that  unless  the  step-size  was  reduced  to  h  <  0.001  the  filter  would  often 
update  Pxx  and  other  matrices  to  non-positive  definite  matrices  (non- valid  updates).  While 
decreasing  the  step-size  did  reduce  the  number  of  non- valid  updates,  the  incidence  of  non- 
valid  updates  was  not  completely  removed.  When  this  positive  definite  property  is  lost, 
the  filter  quickly  diverges. 
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Alternative  implementations,  that  were  not  tried,  may  have  reduced  the  incidence  of  non- 
valid  updates  these  include:  reformulation  of  the  update  equations  to  maintain  the  positive 
definiteness  property,  and  the  use  of  other  numeric  integration  techniques. 

Instead,  the  implemented  code  checks  for  non- valid  updates  of  the  matrices  and  attempts 
to  correct  the  updates  and  keep  the  filter  stable.  In  simulation  studies  in  Section  6, 
simulations  involving  situations  where  the  PME  had  extensive  numeric  problems  have 
been  excluded.  Hence,  the  reader  must  consider  the  results  as  best  case  results,  indicative 
of  the  sort  of  performance  that  may  be  possible  for  a  hybrid  system  filtering  approach, 
rather  than  an  endorsement  of  the  PME  filter. 

The  numeric  stability  problems  on  the  PME  were  found  to  be  significant  enough  to  exclude 
use  of  the  filter  in  a  practical  environment. 

Target  Motion  Assumptions 

The  target  is  assumed  to  exhibit  piece-wise  perfect  linear  and  perfect  turning  motion.  The 
turn  rate  of  the  target  was  assumed  to  be  perfectly  known  by  the  filter,  but  the  timing  of 
turning  events  is  not. 

These  target  motion  assumptions  as  quite  unrealistic.  In  real  problems  the  turning  rates 
will  neither  be  known  nor  be  constant.  An  extensive  investigation  of  situations  involving 
model  mismatch  is  required  before  any  hybrid  system  filtering  approach  could  be  applied  to 
a  real  problem.  Such  examinations  have  been  performed  for  the  IMM  filter  [9].  However,  in 
this  study,  the  computation  effort  required  to  implement  the  PME  filter  was  too  significant 
for  us  to  test  a  large  set  of  engagements  and  hence  performance  during  mismatch  was  not 
examined. 

Inteceptor  Motion  Assumptions 

A  simple  augmented  proportional  navigation  guidance  law  without  either  guidance  or 
autopilot  lags  was  used  to  evaluate  the  filters  as  part  of  a  guidance  loop. 

Matlab™  Implementation 

Details  of  the  Matlab™  implementation  of  the  filters  are  presented  in  Appendix  A 
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6  Simulation  Studies 

In  this  section  we  examine  the  three  presented  filters  in  an  engagement  against  a  target 
performing  a  series  of  manoeuvres.  We  first  consider  the  error  performance  of  the  filters 
(under  no  guidance)  and  then  we  investigate  the  miss  distance  performance  when  these 
filters  are  used  in  a  missile  guidance  loop. 


6*1  Stand  Alone  Filter  Performance 

To  examine  the  performance  of  the  filters,  the  EKF,  IMM  and  PME  were  used  to  filter  data 
in  a  total  of  1249  simulations.  Details  of  the  parameter  values  used  in  these  simulations 
are  provided  in  Appendix  A.  The  EKF  with  knowledge  of  the  actual  target  mode  was 
also  simulated.  Because  all  the  filters  in  this  study  axe  EKF  based  the  EKF  filter  with 
knowledge  of  the  actual  target  mode  provides  a  lower  bound  on  performance. 

As  discussed  in  the  previous  section,  the  PME  filter  sometimes  has  numeric  difliculty  in 
the  Pxx  update  step.  The  implemented  version  of  the  PME  filter  monitors  the  Pxx  update 
step  and  provides  a  count  of  the  number  of  bad  update  steps.  The  simulation  studies 
suggest  that  when  more  than  5  bad  Pxx  updates  occur  the  performance  of  the  Pxx  can  be 
quite  poor. 

Although,  the  performance  of  the  PME  was  quite  unsatisfactory,  some  selected  results  are 
presented  in  the  tables  below  to  demonstrate  what  performance  gains  might  be  possible  if 
the  numeric  problems  can  be  resolved.  In  Table  1,  the  average  error  performance  of  the 
filters  are  compared  on  the  75  data  sets  for  which  the  PME  had  no  bad  Pxx  updates. 

The  average  mean  error  values  presented  in  the  table  is  the  average  error  over  the  last  100 
points  of  the  selected  simulation  run.  The  variance  of  error  supplied  in  the  table  is  the 
variance  of  individual  simulation  errors  across  the  selected  set.  In  Table  2,  the  average 
error  performance  of  the  filters  are  compared  on  the  244  data  sets  for  which  the  PME  had 
one  or  less  bad  Pxx  updates.  The  other  1005  runs  have  greater  than  one  bad  update. 
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Table  1:  Filter  results  with  no  bad  Pxx  updates 


Filter 

Average  Mean  Error 

Variance  of  Error 

EKF  (truth) 

33.1  m 

82.9  m1 

EKF 

569.0  m 

446.9  m2 

IMM 

118.3  m 

1380.9  m2 

PME 

72.6  m 

332.9  m2 

Table  2:  Filter  results  with  one  or  less  bad  Pxx  updates 


Filter 

Average  Mean  Error 

Variance  of  Error 

EKF  (truth) 

32.6  m 

70.2  m? 

EKF 

574.7  m 

1400.2  m2 

IMM 

116.0  m 

992.7  m 2 

PME 

75.9  m 

1257.3  m2 

6.2  Filters  in  a  Missile  Guidance  Loop 


To  evaluate  the  performance  of  the  PME  and  the  influence  of  having  pose  measurements 
of  the  target  available  in  the  missile  guidance  problem,  the  performance  of  the  PME  with 
various  amounts  of  target  pose  information  was  compared  to  the  performance  of  a  guidance 
loop  using  an  EKF  that  knows  the  present  target  manoeuvre. 

Table  3  shows  the  achieved  miss-distance  in  5  guidance  loop  configurations.  Details  of  the 
parameter  values  used  in  these  simulations  are  provided  in  Appendix  A.  To  examine  how 
the  availability  of  mode  measurements  influence  the  performance,  the  PME  was  examined 
with  four  different  amounts  of  mode  information.  The  scenarios  considered  were: 


•  the  PME  with  mode  measurements  for  the  complete  simulation, 

•  the  PME  with  no  mode  measurements  available, 

•  the  PME  with  mode  measurements  available  when  range  was  less  than  1500  m,  and 

•  the  PME  with  mode  measurements  available  when  range  was  less  than  750  m. 

The  achieved  miss-distances  and  the  count  of  Pxx  bad  updates  are  provided  in  Table  3. 

Due  to  the  large  number  of  Pxx  bad  updates  it  is  not  possible  to  determine  the  influence 
of  mode  measurement  on  guidance  performance. 
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Table 

Filter 

8:  Guidance  Results 

Miss  distance  Pxx  Bad  Updates  | 

EKF(truth) 

0.8514  m 

n/a 

PME(full  mode) 

7.9520  m 

2 

PME(no  mode) 

30.6579  m 

5500 

PME(1500m  mode) 

9.0632  m 

5501 

PME(750m  mode) 

10.6382  m 

5500 

6.3  Simulation  Conclusions 

The  numeric  problems  of  the  PME  axe  a  very  serious  limitation.  Due  to  the  very  long  time 
it  took  to  perform  the  simulations  (2  weeks  for  a  full  simulation  set)  it  was  not  possible 
to  work  towards  a  numerical  stable  implementation  of  the  filter.  However,  by  looking  at 
the  performance  of  the  filter  in  the  sub-set  of  simulations  that  Pxx  was  updated  correctly, 
it  seems  the  filter  may  offer  performance  advantages  over  the  IMM  and  EKF  approaches. 
This  may  also  suggest  that  mode  measurements  can  involve  target  filter  performance. 

The  simulations  on  the  performance  of  the  filter  in  a  missile  guidance  problem  were  not 
conclusive  because  of  the  numeric  problems.  The  study  done  here  was  not  complete  and 
significant  work  is  required  on  the  PME  filter  before  it  could  be  considered  fully  examined. 

7  Conclusion 

This  report  provides  a  description  of  recent  work  examining  the  use  of  hybrid  filters  in¬ 
cluding  the  Polymorphic  Estimator  (PME)  for  the  estimation  of  manoeuvring  targets. 
Details  of  the  filtering  problem  and  standard  approaches  such  as  the  extended  Kalman 
filter  (EKF)  and  the  interacting  multiple  model  (IMM)  filters  were  provided.  The  value 
of  mode  measurement  (or  target  pose  information)  was  examined  and  simulation  studies 
described. 

The  PME,  as  presented  in  the  literature,  was  found  to  have  very  serious  numeric  problems 
that  hindered  a  full  evaluation  of  the  hybrid  filtering  approach.  The  simulation  studies 
suggested  that  the  PME  (and  the  use  of  mode  measurements)  may  offer  performance 
advantages  over  the  well  known  EKF  and  IMM  approaches;  however,  significantly  more 
work  (including  the  development  of  numerically  stable  versions  of  the  filter)  is  required 
before  the  examination  of  the  PME  could  be  considered  complete. 
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Appendix  A  Simulation  Environment 

In  these  appendices  we  describe  the  operation  of  the  filter  in  terms  of  the  notation  used 
to  implemented  the  filter  in  the  MatlabTAf .  The  Matlab™  code  was  used  to  examine  the 
PME,  IMM  and  EKF  filters.  The  filters  were  examined  in  two  problems:  As  a  stand  alone 
filter  (representing  measurements  from  the  origin),  and  as  part  of  a  missile  guidance  loop. 
Below  we  provide  a  discussion  of  the  implementation  of  each  simulation  separately. 

A.l  Path  generation:  Overview 

A  trajectory  of  a  missile  is  simulated  in  Mat.lab™.  At  each  time  instant,  the  kinematic 
state  of  the  missile  is  represented  by  a  state  vector,  [x;y]vx;vy],  (x,y  position  and  x,y 
velocity)  which  are  stored  as  columns  in  the  array,  Xk. 

The  missile  in  the  simulation  performs  a  series  of  manoeuvres,  including  straight  sections 
and  turning  left  or  turning  right  segments.  In  the  code  presented  above  the  missile  target 
performs  5  manoeuvres  on  its  flight  toward  the  origin  : 

•  turn  left  (5sec), 

•  straight  (6sec), 

•  turn  right  (5sec), 

•  straight  (3sec)  and 

•  turn  left  (5sec). 

The  target  begins  at  an  initial  state,  xs  (with  units  m  and  m/s)  and  all  turns  are  with 
the  same,  ±phi  turn  rates  (units  rad/s). 

Simulated  measurements  of  the  target  from  the  coordinate  origin  are  also  generated,  in¬ 
cluding: 

•  A  simulated  radar  producing  noisy  measurements  of  range  and  bearing. 

•  A  simulated  imager  (imagerl.m)  that  measures  the  target  orientation  under  noise. 
Target  orientation  is  specified  to  a  30°  sector  (bin)  numbered  1-12  from  North. 
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The  radar  measurements  and  orientation  are  combined  and  produce  a  vector  of  measure¬ 
ments  and  are  stored  in  Ymk  ([range;bearing;bin]) ,  The  mode  (or  current  manoeuvre)  of 
the  target  is  recorded  as  modev . 

A*2  Path  Generation:  Sub-functions 

In  the  implementation,  the  whole  trajectory  of  the  missile  is  simulated  by  combining 
simulated  trajectories  for  each  distinct  manoeuvre.  For  each  manoeuvre  the  function 
asm.m  generates  a  discrete  series  of  states  ( xkl )  with  sampling  time  (ft),  typically  0.001s, 
Internally,  the  complete  path  state  record  is  stored  as  the  vector  Sk,  Measurements  are 
available  with  sampling  time  (fcm),  typically  0.1  s,  in  the  vector  Ymk. 

The  initial  state  ( xkm )  is  evolved  at  each  time  instant  k  using  the  following  state  equation: 

xkl  —  Ah  x  xkm  +  G  x  w 

where  xkm  is  the  state  vector  for  the  previous  sample  time;  Ah  is  a  4  x  4  matrix  selected 
to  produce  a  path  which  is  straight  or  turning  left  or  right  at  typically  0.2rad/s;  G  is  the 
covariance  matrix  of  the  noise  components;  and  w  is  the  process  noise,  set  to  zero  for  all 
simulations. 

At  each  measurement  instant,  i,  a  rectangular  -  polar  coordinate  conversion  is  performed 
giving  the  target  range 

Ck(l,i)  —  x  ns))2  4-  (Sk(2,i  x  ns))2; 

and  bearing 

Ck{ 2,  i)  =  aian(Sk( %  %  x  ns),  Sk(  1,  i  x  ns)); 

where  ns  -  hmjh  must  be  an  integer  and  i  is  the  measurement  interval  (hm)  index.  To 
these  range  and  bearing  values  measurements  noise  is  added: 

Vk(Yi)  =  \frvar  x  randn;  (Range  measurement  noise) 

and 

Vk( 2,  i)  ~  yjtvar  x  randn ;  (Bearing  measurement  noise) 

Here  war  is  the  range  variance,  tvar  is  the  bearing  variance  and  randn  is  Matlab’s 
normal  distributed  random  number  generation.  The  radar  measurements  are  recorded  in 
Ymk(  1  :  2 ,i)  =  Ck(:,i)  +  Vk(:,i);, 
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For  each  measurement,  i,  imagerl.m  works  out  the  orientation  of  the  target  (0)  and  assigns 
a  bin  number  (bin).  Imager  errors  are  introduced  randomly  to  simulate  the  following  : 

•  Total  occlusion  (UDE  -  5%  chance)  -  selects  any  bin  number  at  random. 

•  Nearest  neighbour  (NNE  -  5%  chance)  -  assigns  one  of  the  bins  either  side  of  the 
correct  one. 

The  output  bin  is  recorded  as  Ymfc(3,i). 

At  the  completion  of  each  manoeuvre  the  state  record  ( Xk ),  observations  (Ymfc),  mode 
record  (modev  1)  are  stored.  The  final  state  (xkm)  of  each  manoeuvre  segment  becomes  the 
initial  state  of  the  next  manoeuvre.  The  successive  Xfc,  Ymk  and  modevl  are  concatenated 
to  produce  a  complete  state,  measurement  and  mode  record  for  the  path. 
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Appendix  B  Filtering  Evaluation  Code 

The  following  subsections  describe  the  various  routines  used  to  compare  the  EKF,  IMM 
and  PME  filter  on  the  target  filtering  problem. 

imtcommon.m  Contains  initialization  and  definitions  for  a  number  of  parameters  and 
matrices  which  are  common  to  all  the  filters.  These  include  sample  time  (h),  measurement 
time  ( hm ),  turn  rate  (phi )  and  state  transition  matrices  (Ac,  Ad).  The  following  vectors 
and  matrices  are  also  initialised:  initial  target  state  (xs), initial  target  state  estimate  ( xhs ), 
measurement  and  process  noise  covariance  matrices  ( Rn ,  Qn),  and  the  prior  state  error 
covariance  matrices  ( Pn ). 

The  quantity  Ac  is  (4  x  4  x  3)  3d-matrix  where  the  third  index  corresponds  to  the  selection 
of  a  2d-matrix  representing  the  continuous-time  state  transition  matrices  corresponding 
to  turning  left,  straight  flight,  and  turning  right  respectively.  The  quantity  Ad  is  similarly 
used  to  represent  three  discrete-time  state  transition  matrices. 

The  measurement  and  process  noise  covariance  matrices  (Rn,  Qn)  are  individually  stored 
for  EKF,  IMM  and  PME  filters.  The  prior  state  error  covariance  matrix,  Pn,  is  initialised 
with  errors  of  100  m  in  position  and  10  m/s  in  velocity  for  all  runs  (unless  a  larger  state 
error  is  input  via  xhs). 

EKF  filters  The  EKF  algorithm,  as  described  in  (3.3)  was  implemented  with  notation 
changes  to  allow  for  the  restricted  text  format  options  in  the  Matlab™  code  editor.  These 
changes  are: 


Maths  Terms 

Matlab™  Term 

Notes 

%k\k- 1 

xh 

xhm 

xhm  is  the  output  of  one  iteration  (%|&) 
and  also  the  input  of  the  loop(%_4jk_i). 

Ahm 

^*k\k— 1 

Pm 

■^k\ky  'ffc-lffe-l 

Pm 

Note:  Pm  is  both  the  output  of  one  iteration  (P^jfc) 
and  the  input  of  the  next  loop  ( Pk-i\k-i )* 
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Ck 

C 

Kk 

Kk 

Vk 

Ymk 

Ck&k\k-1 

yhat 

Rk 

R 

Qk 

0 

we  assume  there  is  no  process  noise 


At  each  measurement  instant  the  appropriate  state  transition  matrix  Ad(:, :,i)  is  selected 
(assuming  the  EKF  has  access  to  the  true  target  manoeuvre  state.  The  variable  prf  is 
used  to  switch  between  two  measurement  situations.  If  prf  =  10  then  the  measurement 
update  occurs  at  hm .  If  prf  =  any  other  value  then  the  measurement  update  occurs  at 
10  x  hm.  The  variable  xhatv  is  the  complete  record  of  the  estimated  state  output  of  the 
filter,  at  the  rate  dictated  by  hm. 

For  the  ekfc.m  version,  no  knowledge  of  the  target  manoeuvre  is  assumed,  so  the  state 
transition  matrix  corresponding  to  straight  flight  Ad(:, :,  2)  is  used. 

IMM  filter  The  IMM  filter,  as  described  in  (4.1)-(4.4),  was  implemented  with  notation 
changes  to  allow  for  the  restricted  text  format  options  in  the  Mat  lab™.  These  changes 
are: 


Maths  Terms 

Matlab™  Term 

Notes 

M«| -  1) 

mu2 

Equation  (4.1) 

rt-ilk-i 

xoj 

Equation  (4.1) 

pjo 

*k—l\k—l 

Poj 

Equation  (4.1) 

4\k-l 

xjh 

Equation  (4.2) 

Pk\k-1 

pj 

Equation  (4.2) 

4 

Kk 

Equation  (4.2) 

xhjk 

Equation  (4.2) 

Pk\k 

Pjk 

Equation  (4.2) 

4 

Ljk 

Equation  (4.3) 
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fi 

mu 

Equation  (4.3) 

%k\k 

xkk 

Equation  (4.4) 

Pk\k 

Pkk 

Equation  (4.4) 

It  is  assumed  in  the  initial  definition  of  the  mode  probability  (mu)  that  the  the  target  is 
equally  likely  to  initially  be  in  one  of  the  three  possible  modes  (straight  path,  turning  left 
or  right).  It  is  also  assumed  that  the  target  will  continue  in  the  same  mode  (manoeuvre) 
or  the  target  can  switch,  with  equal  probability,  to  any  adjacent  manoeuvre  mode  (ie.  the 
model  does  not  allow  for  direct  transition  between  left  and  right  turn  modes).  This  mode 
switch  is  random  and  described  by  a  probability  transition  matrix  introduced  in  the  PME 
description  below. 

PME  filter  The  PME  filter,  as  described  in  (4.5)  -(4.8),  was  implemented  with  notation 
changes  to  allow  for  the  restricted  text  format  options  in  the  Mat  lab™.  These  changes 
are: 


Maths  Terms 

Matlab™  Term 

fhm 

mut 

fh 

dx  x  deltaT 

It 

dPxf  x  deltaT 

A  p 

dPxx  x  deltaT 

d  p 

dPxxfm  x  deltaT 

P/ifj, 

Pff 

Pxfl 

Pxf 

p 

Pxffi 

P{xp^}xilm 

Pxfixfm 

P^Xflm 

Pfixfm 

PxXpL*  jlm 

Pxxfifm 

y f 

ybm 

Notes 

Time  update  equations 


Measurement  update  equations 


23 


DSTO-TN-0488 


4 

fh 

xh 

Ax 

dxh 

P r/i 

Pxf 

P rx 

Pxx 

lx 

Kk 

This  filter  is  based 

on  the 

Range/bearing  measurements  update  equations 


Boyd  and  Sworder’s  version  of  the  Matlab™  code  is  incorporated  using  wrapping  code 
to  translate  between  the  different  variable  names  used  in  their  implementation. 


In  defining  a  hybrid  system  description  appropriate  for  the  PME  filter  several  design 
assumptions  were  made.  It  was  assumed  that  the  target  was  able  to  perform  3  types  of 
manoeuvre  (straight  flight,  left  turn,  right  turn).  The  probability  of  transitioning  between 
manoeuvre  modes  is  described  by  the  (3  x  3)  transition  probability  matrix  Aa  .  The 
elements  of  Aa  are  such  that  A(i,j)  represents  the  probability  of  changing  from  mode  j 
to  mode  i. 


At  each  time  instant  the  target  orientation  is  one  of  twelve  values  (corresponding  to  the 
orientation  bins).  The  probability  of  changing  orientation  bin,  when  that  target  is  in  a 
particular  manoeuvre  mode,  is  described  by  3  (12  x  12)  matrices,  As,  Al  and  Ar. 

It  was  decided  to  use  a  composite  discrete  state  vector,  /h,  to  represent  both  the  present 
manoeuvre  mode  and  orientation  of  the  target.  This  vector  represents  all  combinations 
of  the  12  orientations  and  3  manoeuvre  modes  so  fh  is  a  (36  x  1)  vector.  The  transition 
probability  matrix,  A,  required  by  the  PME  is  hence  a  (36  x  36)  matrix  made  up  from 
Aa,  As,  Al  and  Ar. 

The  PME  filter  also  requires  an  observation  mapping  matrix,  P6,  which  describes  the 
probability  of  receiving  ybm  given  that  fh  has  a  particular  value.  Because  the  orientation 
observation  is  independent  of  the  target  manoeuvre,  the  mapping  matrix,  P6,  contains 
three  repeated  blocks,  Po. 
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Appendix  C  Evaluation  of  Guidance 
Performance  Code 

To  evaluate  tlie  filter  as  part  of  a  guidance  loop  we  simulated  a  missile  guidance  problem. 
In  the  missile  guidance  loop,  the  missile  observes  the  target  through  a  range  and  bearing 
measurements.  The  missile  uses  these  measurements  to  estimate  the  state  of  the  target 
and  then  uses  these  state  estimates  to  guide  towards  the  target. 

The  guidance  algorithm  used  in  these  guidance  loop  simulations  is  the  augment  propor¬ 
tional  navigation  algorithm: 

ut  =  -3(FCA  +  id).  (Cl) 

where  d  is  the  target  acceleration  perpendicular  to  the  line-of-sight  line  connecting  the 
missile  and  the  target;  Vc  is  the  closing  velocity;  A  is  the  line-of-sight  rate  seen  by  the 
missile;  and  Ut  is  the  control  action. 

In  previous  sections  we  have  introduced  all  the  information  necessary  to  construct  the 
state  estimates  so  this  is  not  repeated  here.  The  control  algorithm  described  above  was 
implemented  with  notation  changes  to  allow  for  the  restricted  text  format  options  in  the 
Mat  lab™  code.  These  changes  are: 


Maths  Terms 


Matlab™  Term 
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Appendix  D  Matlab™  Syntax  Headers 


In  this  section  we  introduce  the  headers  of  the  Matlab™  code  used  examine  various  filters. 


Stand  Alone  Filters  Figures  D1  and  D2  demonstrate  the  interconnection  of  Matlab™ 
files  used  to  examine  the  filters  in  a  stand  alone  implementation. 


pathlsrsl.m 


i 

i 

Xk  ! 

asm 

_ yk_ 

i 

i 

— i — 

imagerl 

mode 

i 

i 

— i — 1 

i 

i 

i 

Target 
State  and 
Measurements 


Figure  Dl:  Matlab™  files  for  data  generation  in  stand  alone  filtering  simulation 


Path  Generation 


•  pathlsrsl 

•  asm 

•  imager  1 


'/.Syntax: 

'/.Description 

X 

'/.Inputs 
X Out puts 

X 


l 

% 

I 

% 


:  function  [Xk, Ymk, mode v, xho, vk, phi, t,hm,h]=pat hi srsl( dummy) ; 

:  Generates  a  target  path  in  2-d  space (x-y) and  encompasses  6  manoeuvres  - 
left  tum(5sec)  .straight (6sec)  .right  turn (5sec)  .straight (3sec)  .left  turn(Bsec) 

:  nil  (all  necessary  variables  are  set  by  initcommon) 

:  Xk-path  state  vector  (4x(T/h))  (T=total  path  time  (sec),  h=sample  time  (sec)) 
Xk(l,:)«x  position,  Xk(2,:)=  y  position, Xk(3, :)*=x  velocity,  Xk(4,:)=y  velocity 
:  Ymk-radar  measurements 

Ymk(l, :)=range(m) ,  Ymk(2,:)=  bearing (deg.) 

Ymk(3, :)*target  direction  of  travel (bin  number,  of  30deg  opening) 

:  modev-record  of  the  target  mode (straight  path  or  turning  left/right) at 
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Target 

filter 

estimates 

xhatv 


Figure  D2;  MatlabT ^  files  for  filters  in  stand  alone  filtering  simulation 


% 

%  : 

%  z 

% 

% 

%  t 

%  z 

% 

%  z 

XCalls  : 

% Authors  : 

^Modifications  ; 
*/iBugs  ; 

^Assumptions  ; 


each  sample  point 

xho*“  estimate  of  the  state  vector  at  the  beginning  of  the  run 
vk-  process  and  measurement  noise  variance  vector 

vk(l)*  range  measurement  noise  variance, vk( 2) ^bearing  measurement  noise  variance 
vk(3)=x  component  process  noise  variance,vk(4)=  y  component  process  noise  variance 
phi*  turn  rate  (rad/sec)  -  -Mre  =left  turn,-ve  =  right  turn  of  last  manoeuvre. 
t=  Total  flight  time 
hm=measurement  interval  (sec) 
h*sample  interval  (sec) 
asm.m,  init common. m 
Jason  Ford,  Peter  Hunter 
final  version  14-11-01 
nil 

Target  has  constant  velocity  segments  and  constant  turning  rate  segments. 


/(Syntax: 

/(Description 

% 

/(Inputs 

% 

% 

% 

% 

% 

% 

% 

y.  Out  puts 
% 


z  function  EXk,¥mk,xkffi,nm,ffiodevl]=asm<phi,xkm,vk,t,hm,h) ; 

:  Generates  a  target  path  in  2-d  space (x-y) encompassing  1  manoeuvre  - 
left  turn  (phi*+0,2  rad/sec) , straight (phi=0  rad/sec)  or  right  tum(phi=-0.2  rad/sec) 
:  phi*  turn  rate  (rad/sec)  -  +ve  =left  turn,-ve  =  right  turn 
:  xkm  =  initial  state  vector  for  this  manoeuvre  iteration. 

:  vk-  process  and  measurement  noise  variance  vector 

vk(l)=  range  measurement  noise  variance ,  vk  (2) =bearing  measurement  noise  variance 
vk(3)*x  component  process  noise  variance , vk (4) =  y  component  process  noise  variance 
:  t=  manoeuvre  time 
:  hm=measurement  interval  (sec) 

:  h*sample  interval  (sec) 

:  Xk-path  state  vector  (4x(T/h))  (T*total  path  time  (sec),  h*sample  time  (sec)) 
Xk(l,:)=x  position,  Xk(2,:)=  y  posit ion, Xk (3, ;)=x  velocity,  Xk(4,:)*y  velocity 
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•/. 

y. 

X 

'/. 

X 

X 

X 


X Calls 

/(Authors 

*/, Modifications 

'/.Bugs 

’/.Assumptions 


Ymk-radar  measurements 

Ymk(l,  :)=range(m) ,  Ymk(2,:)*  bearing(deg.) 

Ymk(3, :)=target  direction  of  travel (bin  number,  of  30deg  opening) 
xkm  -  (4x1)  initial  state  vector  for  the  manoeuvre 
nm  -  the  number  of  measurements  in  the  manoeuvre 

modevl-record  of  the  target  mode (straight  path  or  turning  left/right) at 

each  sample  point 

imagerl.m 

Jason  Ford,  Peter  Hunter 
final  version  15-11-01 
nil 
nil 


'/.Syntax : 
y, Description 

X 

X 

X 

X 

X 

X Inputs 

X 

X 

X 


X 

X 

/(Outputs 

'/.Calls 

^Authors 

y, Modifications 

'/.Bugs 

/(As sumptions 

X 


[bin]  *=imagerl  ( Ymk , Vy , Vx , Y , X , i) ; 

The  target  direction  is  located  in  one  of  twelve  30  deg. intervals (labelled 

cw  from  North  as  {1,2,3, .. ,12}) 

Randomly  occurring  allocation  errors  are  introduced: 

-NNE  -  nearest  neighbour  error-  placed  into  an  adjacent  bin  due  to  noise  etc. 
-UDE  -  uniformly  distributed  error  -  occlusion  of  target  makes  any  bin  number 

a  possibility 

Ymk-not  used 

Vy  -  target  velocity  -  y  direction 
Vx  -  target  velocity  -  x  direction 
Y  -  target  location  -  y  direction 
X  -  target  location  -  x  direction 
i  -  not  used 

bin  -  the  direction  of  target  travel  as  classified  {1:12}  (noisy  result) 
nil 

Jason  Ford,  Peter  Hunter 
final  version  14-11-01 
nil 

Target  can  be  represented  by  1  of  twelve  attitude  bins. 

Assumptions  on  the  nature  of  the  error. 


Filters 


•  pme 

•  imm 

•  ekf 


'/.Syntax : 
'/.Description 

X 

X 

'/.Inputs 

'/. 

X 

X 


/(Output 

X 


:  [xhatp,Pxxv,normrecord,modeestv,binestv]=pme(Ymk,Xk) ; 

:  implements  a  Polymorphic  Estimator  algorithm  to  produce  an  estimate 
of  position  and  velocity  of  the  target  derived  from  noisy  range 
and  bearing  measurements  from  a  simulated  tracking  radar. 

:  Ymk-radar  measurements  from  the  simulation  origin. 

Ymk  ( 1 ,:)  grange  (m) ,  Ymk(2,:)B  bearing  (deg.) 

Ymk(3, :)*=target  direction  of  travel(bin  number,  of  30deg  opening) 

:  Xk-target  state  (x,y,x-velocity,y-velocity)  and  only  used  for  the  initial  Pxxv 
:  xhatp  (4xT)  vector  of  position  estimate  (x»y,  x-velocity,  y-velocity) 

:  Pxxv (4xT) -vector  of  Pxx  diag  elements  (used  by  plotdata  routine). 
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/ 

% 


% 

‘/.Calls 

'/♦Authors 

/Modifications 

/.Bugs 


% 

% 

Z 


/Assumptions 


:  normrecord-  debuging  info. 

*  modeestv(lxT) -vector  of  manu.  mode  estimates. 

;  binestv(lxT) -vector  of  orientation  mode  estimates. 

:  initpme 

:  Jason  Ford,  Peter  Hunter 
:  final  version  22-11-01 

:  The  time  update  step  of  the  algorithm  is  not  particularly  stable. 

Small  time  step  sizes  (ie.  h)  seem  to  improve  the  preformance  but 
bad  updates  of  Pxx  vector  are  possible  (if  this  happens  the  code  resets 
some  quantities  in  the  filter  to  try  to  stop  the  filter  diverging) . 

:  numeric  intergration  of  continuous-time  equations. 


/Syntax: 

'/Description 

Z 

% 

/ 

Z 

/Inputs 

Z 

% 

% 

/Outputs 

% 

% 

Z 

% 

% 

/Calls 

/Authors 

/Modifications 

/Bugs 

/Assumptions 


:  (xhat i ,  Pkkv  ,muv ,  xhl ,  xh2 ,  xh3]  =imm(Ymk) ; 

;  Interacting  Multiple  Mode  filter  (IMM)-  is  a  filter  function  used  to  estimate 
the  location  of  a  manoeuvring  target.lt  uses  measurement  data  as  veil  as  manoeuvre 
type  (straight  or  turning  path)  to  estimate  the  path  of  the  target. 

It  uses  the  mode  estimation  to  switch  between  individual  EKF  algorithms  to  better 

track  the  target. 

:  Ymk-radar  measurements  from  the  simulation  origin. 

Ymk(l , : )  =range (m) ,  Ymk<2,;)=  bearingCdeg.) 

Ymk(3, :)=target  direction  of  travel (bin  number,  of  30deg  opening) 

.modev  (IxT/h)  vector  which  logs  the  actual  manoeuvre  at  each  sample  point 
:xhati  (4xT)  vector  of  estimates  (x,y,  x-velocity,  y-velocity) 

;Pkkv  combined  state  covariance  matrix 
:muv  mode  probability  for  each  modeal  filter 
:xhl  state  estimate  of  filter  1 
:xh2  state  estimate  of  filter  2 
:xh3  state  estimate  of  filter  3 
:  initcommon 

:  Jason  Ford,  Peter  Hunter 
;  final  version  20-1 i-01 
:  nil 

:  Based  on  hybrid  system  mode  of  dynamics. 


/Syntax: 

/Description 

Z 

Z 

Z 

/Inputs 

z 

z 

•/ 

z 

/Output 

/Calls 

/Authors 

/Modifications 

/Bugs 

/Assumptions 


function  [xhat v] =ekf (Ymk ,modev , prf ) ; 

implements  an  Extended  Kalman  Filter  to  produce  an  estimate  of  position  and 
velocity  of  the  target  derived  from  noisy  range  and  bearing  measurements  from 
a  simulated  tracking  radar.  The  filter  has  the  target  manoeuvre  at  each 
sample  point  (straight  path  or  turning  left  or  right  at  0.2  rad/sec) 

Ymk-radar  measurements  from  simulation  origin 

Ymk(l* :)=range(m) ,  Ymk(2,:)=  bearingCdeg.) 

Ymk (3, :)=target  direction  of  travel (bin  number,  of  30deg  opening) 
modev  (IxT/h)  vector  which  logs  the  actual  manoeuvre  at  each  sample  point 
prf  -  measurement  interval  (sec)  (pulse  repetition  frequency  of  the  radar  system) 
xhatv  (4xT)  vector  of  position  estimate  (x,y,  x  velocity,  y  velocity) 
initcommon 

Jason  Ford,  Peter  Hunter 
final  version  20-11-01 
nil 
nil 


Initialisation  files 

•  initcommon 
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•  initpme 


‘/♦Syntax: 
•/♦Description 
X Inputs 
•/♦Outputs 

7. 

% 

7. 

7. 

7. 

7. 

•/.Calls 
•/♦Authors 
•/.Modi  f  i  c  at  i  ons 
•/.Bugs 

•/.Assumptions 


:init common 

:  contains  the  initialization  values  for  all  filters 
:  nil 

:  there  are  no  explicit  outputs  but  the  following  variables  are  initialised: 

-  h  (sample  interval (s) ), phi (turn  rate (rad/sec) ) ,hm 

-  xs,xhs 

-  Ac, Ad,Rn,Qn,Pn 

-  for  EKF  :  EKFQ,EKFR,EKFP 

-  for  IMM  :  IMMQ , IMMR , IMMP 

-  for  PME  :  Rx,Rxx,Pxx 
:  nil 

:  Jason  Ford,  Peter  Hunter 
:  final  version  14-11-01 
:  nil 
:  nil 


•/.Syntax: 

•/.Description 

•/.Inputs 

•/.Outputs 

7. 

7. 

7. 

7. 

•/.Calls 

•/.Authors 

•/.Modifications 

‘/.Bugs 

‘/.Assumptions 


initpme 

Contains  constants  and  initialization  data  for  the  pme.m  program 
nil 

there  are  no  explicit  outputs  but  the  following  variables  are  initialised: 

-  statenum,modenum,binnum,measbin, bigmode 

-  hm,minfh 

-  Aa,As,Al,Ar»Ai 

-  Po,Q,Rx,Rxx. 
ini t common 

Jason  Ford,  Peter  Hunter 
final  version  27-11-01 
nil 
nil 


Guidance  Evaluation  Figure  D3  describes  the  structure  of  the  Mat  lab™  code  used 
to  simulate  a  missile  guidance  loop  using  an  extended  Kalman  filter.  The  Matlab™  code 
for  incorporating  the  PME  filter  is  similar. 


Matlab  Files 

•  comparefilters 

•  guidanceekf 

•  guidamcepme 

•  ekfscon 

•  pmelloop 
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initekf.m 

first.m 


Manoeuvre  Defn 


Guidance  Law 


State  step  forward 


Observations 


Generated  |  imagerl  ,m 


% 


mode 


ekfscon 


Update  Output 
Matrices 


End  Engagement 


xhatv 
xhatvp 
Xk 
Xv 
Ymk 

**►  minv 


Figure  D3 :  Matlab™  files  for  guidance  loop  simulation  using  the  extended  Kalman  filter . 


•  imagerl 

•  first 

•  initekf 


•  initpme 


Kcomparef liters 
*/,  Syntax:  : 

%Description  : 
^Inputs  ; 

*/,  Out  puts  : 

%  t 

%  : 

% 

% 


% 


% 

'/.Calls 


comparefilters 

runs  a  number  of  control/filter  algorithms  to  compare  the  minimum  miss  distances  (minv) 
nil 

MISSmin  -  vector  of  the  minimum  miss  distances  for  each  filter 
MeanMISSminv  -  for  multiple  runs  this  is  the  average  result  for  each  filter 
additional  outputs  with  prefixes  :EKF. . .  (c)  tPMEM0B, , .  ,PMEKM. . ,  .PHESM1500. 

PMESKT50. . . 

-  ekf  filter  ;Pm,Xk,xhatv>Xv,€latv,Tlatv 

-  pme  f liter  : Pxxv  *  Xk , xhatv , pxxcount  *pxxvcount , badextraps , Xv , Clat v , Tlatv 
guidanceekf ,  guidancepme,  first 
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'/.Authors 

'/.Modifications 

'/. 

X 

X 

'/. 

X 

% 

% 

% 

•/.Bugs 

'/.Assumptions 


Jason  Ford,  Peter  Hunter 


16/2/01 

21/2/01 

28/2/01 

12/7/01 

18/7/01 

2/10/01 

3/10/01  i 
27/11/01: 
:  nil 
:  nil 


save  workspace  to  comparef4a.mat 

save  workspace  to  comparef4b.mat 

save  workspace  to  comparef4c .mat  h=0.0001  s=3 

save  workspace  to  comparef4f .mat  h=0.001,  ekf,  ekfc,  pmemod 

save  workspace  to  comparef 4g . mat  h=0.001  s=10  all  except  ekfc 

save  workspace  to  comparef 4a2. mat  h=0.0001,  s«l,  run=10,  only  EKF  ft 

PMEMOD 

save  workspace  to  comparef 4a3. mat  as  for  comparef 4a2  except  all  filters 
final  version  ,  renamed  comparef liters 


*/,  Syntax: 

*/,De  script  ion 
*/. 

X 

'/,  Global  s 
'/.Inputs 

X 

*/,  Out  puts 
'/. 

*/. 

'/. 

t 

'/. 

% 

'/. 

'/. 

'/. 

1 

*/. 

'/. 

X 

*/, Calls 

X 

X 


*/. Authors 

/.Modifications 

'/.Bugs 

*/,  Assumptions 


[Pmt ,  Xv ,  Xk ,  Ski ,  Ymk ,  mode  v ,  xhat  v ,  xhat  vp ,  Range  v ,  Tlat  v ,  minv ,  Clatv]  *guidanceekf  (time  steps, tiO); 
A  simulated  engagement  between  an  incoming,  manoeuvring  target  and  a  radar  guided 
interceptor.  The  radar  range  and  bearing  measurements  are  filtered  by  an  EKF,  with 
the  whole  operating  vithiu  an  APN  guidance  loop. 

B,  Q,  R,  Bt; 

timesteps  -  the  number  of  time  sample  points  *  realtime/h. 
tiO  -  interceptor  start  heading  (rad) 

Pmt  -  diagonal  elements  of  the  covariance  matrix.  It  is  a  diagnostic  output  to 
enable  construction  of  the  1  sigma  circles  at  each  measurement  point 
of  the  filter  output. 

Xv  -  relative  state  between  interceptor  and  target 

Xk  -  target  state 

Ski  -  actual  target  state 

Ymk  -  measurement  sequence  (radar  range  and  bearing  measurements) 

modev  -  the  target  mode  vector 

xhatv  -  estimate  of  target  state 

xhatvp  -  one-step  ahead  estimate  of  target  state 

Rangev  -  target- interceptor  separation (used  for  testing  guidance  loop) 

Tlatv  -  perpendicular  component  of  the  target  acceleration (testing  variable) 

minv  -  minimum  miss  distance  between  the  target  and  interceptor 

Clatv  -  Vector  of  guidance  commands 

initekf 

ekf s con. m 

imager l.m 

Jason  Ford,  Peter  Hunter 
final  version  27-11-01 
nil 

interceptor  position  is  prefectly  known...  in  real  case  would  need  INS  solution 


'/.Syntax : 

X 

'/.Description 

X 

X 

*/.  Globals 

X 

X 

X 

X 

X 

X 

/,  Inputs 
*/. 
x 


:  [Xv , Xk , Pxxv , Ski , Ymk , modev , xhatv , xhatvp , Rangev , T1 at v , minv , pxxc ount , pxxmc ount , 
badext raps , Clatv] ^guidancepme (timesteps ,ti0 ,moderangemax) ; 

:  A  simulated  engagement  between  an  incoming,  manoeuvring  target  and  a  radar 
guided  interceptor.  The  radar  range  and  bearing  measurements  are  filtered 
by  a  PME,  with  the  whole  operating  within  an  APN  guidance  loop. 

:  B,  Q,  R,  Bt 

:  TIMEUPDATE  MODEUPDATE  BASEUPDATE  RADAR. ALL-TIME 
:  xh  fh  Pxx  Pxf  Pxxfi  xhm  fhm  Pxxm  Pxfm  Pxxfim  modest 
:  A  Rx  Rxx  V  Q  Pb 
:  pxxcount  badextraps  pxxmcount 
:  bigmode  E 
:  h  P.preset 

:  timesteps  -  the  number  of  time  sample  points  *  realtime/h. 

:  tiO  -  interceptor  start  heading  (rad) 

;  moderangemax  -  the  maximum  distance  (m)  over  which  the  imager  provides  information 
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'/♦Outputs 

% 

/ 

Z 

Z 

% 

% 

% 

Z 

% 

% 

% 

% 

% 

% 

% 

/Calls 

Z 

% 

/Authors 

/Modifications 

‘/Bugs 

‘/Assumptions 


Xv  -  relative  state  between  interceptor  and  target 
Ik  -  target  state 

Pxxv-  diagonal  elements  of  the  covariance  matrix.  Used  to  construction 
one  sigma  circles  at  each  measurement  point  of  the  filter  output. 

Ski  -  actual  target  state 

Ymk  -  measurement  sequence  (radar  range  and  bearing  measurements) 
modev  -  the  target  mode  vector 
xhatv  -  estimate  of  target  position 
xhatvp  -  Debug  variable 

Rangev  -  target-interceptor  separation(used  for  testing  guidance  loop) 

Tlatv  -  perpendicular  component  of  the  target  acceleration  (testing  variable) 

minv  -  minimum  miss  distance  between  the  target  and  interceptor 

pxxcount  -  Error  count 

pxxmcount  -  Error  count 

badextraps  -  Error  count 

Clatv  -  Vector  of  guidance  commands 

initpme 

pmelloop.m 

imager 1 .m 

Jason  Ford,  Peter  Hunter 
final  version  27-11-01 
nil 

interceptor  position  is  prefect ly  known...  in  real  case  would  need  INS  solution 


*Z  Syntax : 

/Description 

'/Globals 

/Inputs 

% 

% 

% 

% 

% 

% 

% 

% 

% 

/Outputs 

Z 

X 

'/Calls 

/Authors 

/Modifications 

'/Bugs 

/Assumptions 


i  [xhp, xhm,Pm]=ekf scon (Ahm,xhm,Pm, ymk, radarpoints,i,u,xik)  ; 

:  EKF  section  -  one  iteration  only 
:  B,  0,  R,  Bt; 

:  Ahm  -  initial  state  transition  matrix 
:  xhm  -  initial  target  state  (k-l|k-l) 

:  Pm  -  initial  covariance  martix  (k-llk-l) 

%  Ymk-radar  measurements  from  the  simulation  origin. 

Ymk(l,  ;>=range(m) ,  Ymk (2, :)=  bearing(deg.) 

Ymk (3 , ; ) —target  direction  of  travel (bin  number,  of  30deg  opening) 
:  radarpoints  -  helps  determine  if  a  measurement  update  is  required 
;  i  -  time  k 

*  u  -  control  action  of  the  interceptor. 

:  xik  -  initial  interceptor  state  (k-l(k-l) 

•  xhp  -  state  vector  -  after  time  update  (k-l|k) 

;  xhm  -  measuremnet  updated  state  (kjk) 

:  Pm  -  initial  covariance  matrix  for  (k+i) 

:  nil 

:  Jason  Ford,  Peter  Hunter 
:  final  version  27-11-01 
:  nil 
;  nil 


/Syntax; 

/Description 

/Globals 

X 

% 

X 

% 

% 


/Inputs 

X 


/ 


;  [dummy]  -pmelloop  (Ymk  ,k , xik, radarpoints  ,modepoints  ,moderangemax , Range)  ; 
;  is  a  single  time  iteration  of  the  PME  filtering  algorithm 
:  TIMEUPDATE  MODEUPDATE  BASEUPDATE  RADAR_ALL_TIME 
;  xh  fh  Pxx  Pxf  Pxxfi  xhm  fhm  Pxxm  Pxfm  Pxxfim  modest 
;  A  Rx  Rxx  W  Q  Pb 
;  pxxcount  badextraps  pxxmcount 
:  bigmode  E 
:  h  P_preset 

;  Ymk-radar  measurements  from  the  simulation  origin. 

Ymk (1 , : ) -range (m) ,  Ymk(2,;)=  bearing (deg.) 

Ymk (3 , : ) =target  direction  of  travel (bin  number,  of  30deg  opening) 
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X 

% 

% 

X 

'/.Outputs 

'/.Calls 

'/.Authors 

'/.Modifications 

'/.Bugs 

*/, Assumptions 

'/.Syntax : 

'/.Description 

*/. 

X 

% 

X 

X 

'/.Global  s 
'/.Inputs 

X 

% 

% 

X 

X 

'/.Outputs 

'/.Calls 

'/.Authors 

'/.Modifications 

'/.Bugs 

X Assumptions 


*/, Syntax : 
'/.Description 

X 

'/.Globals 
'/.Inputs 
*/, Out  puts 

V. 

X 

X 

X 

X 

X 

'/.Calls 
V, Authors 
'/.Modifications 
'/.Bugs 

'/.Assumptions 


'/.Syntax : 

'/.Description 

XGlobals 

'/.Inputs 

'/.Outputs 

% 


:  radarpoints  -  number  of  radar  measurements 
:  modepoints  -  number  of  mode  observations 

:  moderangemax  -  maximum  range  (m)  over  which  the  image  information  is  available 
:  Range 

:  nil  explicitly  (output  via  global  variables) 

;  nil 

:  Jason  Ford,  Peter  Hunter 
:  final  version  27-11-01 

:  Same  bugs  as  PME  filter.  Numerical  unstable. 

:  Same  as  the  PME  filter. 

:  [bin]=imagerl(Ymk,Vy,Vx,Y,X,i) ; 

:  The  target  direction  is  located  in  one  of  twelve  30  deg. intervals (labelled  cw 

from  North  as  {1,2,3, .. ,12» 

Randomly  occurring  allocation  errors  are  introduced: 

-NNE  -  nearest  neighbour  error-  placed  into  an  adjacent  bin  due  to  noise  etc. 
-DDE  -  uniformly  distributed  error  -  occlusion  of  target  makes  any  bin  number 
a  possibility 

:  nil 

:  Ymk-not  used 

:  Vy  -  target  velocity  -  y  direction 
:  Vx  -  target  velocity  -  x  direction 
:  Y  -  target  location  -  y  direction 
:  X  -  target  location  -  x  direction 
:  i  -  not  used 

:  bin  -  the  direction  of  target  travel  as  classified  {1:12}  (noisy  result) 

:  nil 

:  Jason  Ford,  Peter  Hunter 
:  final  version  27-11-01 
:  nil 
:  nil 


: first 

:  contains  the  initialization  values  for  all  filters  evaluated  in  the  control  loop 
environment 


nil 

nil 

there  are  no  explicit  outputs  but  the  following  variables  are  initialised: 
-  Basic  System  Parameters  :h  ,phi,hm,modeh,radarh 


-  Sample  Numbers/sec 

-  Initial  State  vectors 

-  Noise  variances 

-  Debugging  Switches 

-  Iteration  Parameters 
nil 

Jason  Ford,  Peter  Hunter 
final  version  27-11-01 
nil 
nil 


:  modepoints , radarpoints , contr olpo int s 

:ta,xs,tiO,xis 

: pxvar , py var , rvar , t var 

: TIMEUPD ATE , MODEUPD ATE , BASEUPDATE , RADAR  ALL  TIME 
: s , runnum , seedof f 


:initekf 

:  contains  the  initialization  values  for  EKF  filter 


:  B,  Q,  R 
:  nil 

:  there  are  no  explicit  outputs  but  the  following  variables  are  initialised: 
-phi 
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%  -xhs 

*/•  -Ac ,  Ad ,  Aint ,  Aintd , Bt ,  Bi 

%  -Pn,Qn,Rn 

%  -ekfp,ekfq,ekfr 

'/.Calls  ;  nil 

/(Authors  :  Jason  Ford,  Peter  Hunter 

^Modifications  :  final  version  27-11-01 
‘/.Bugs  :  nil 

% Assumptions  :  nil 


*/,  Syntax: 

/(Description 

%Globals 


/(Inputs 

^Outputs 

% 

% 

% 

% 

% 

% 


% 

% 


% 


% 

KCalls 

^Authors 

/(Modifications 

%Bugs 

/(Assumptions 


:  initpme 

:  contains  the  initialization  values  for  the  guidancepme  program 
:  ail 
;  nil 

*  there  are  no  explicit  outputs  hut  the  following  variables  are  initialised: 
”s  tat  enum  ,modenum ,  binnum ,  measbin ,  bigmode 
-pxxcount 
-f  ,fh,f  hm,minf  h 
-v 

-xhm , xhs  >  xh , xkm , xo , xho , ym , ybm 

-Aa,As,Al,Ar,A,Ah,Ahm,Ac,Ad,Ai,Aint,Aintd 

-Bt,Bi,C,Ck,E,G 

-Pn,Po ,Pbl ,Pb  ,Pxf ,Pxffi ,Pxx , Q , Qn >Rn ,Rx ,Rxx 
-MKL  >Pxx_p,pixx_p  ,Pxx_m ,  pixx_m ,phi_m ,  phi_p ,  Pxp_m, 

-Pxp_m,Pxp_p,x_p , x_m .badextraps 
:  first 

:  Jason  Ford,  Peter  Hunter 
:  final  version  27-11-01 
:  nil 
:  nil 


35 


DSTO-TN-0488 


Appendix  E  Simulation  Parameter  Values 

E.l  Stand-Alone  Filter  Evaluation 


The  following  parameter  values  were  used  in  the  stand-alone  filter  evaluation: 


h  =  0.001s 
hm  —  0.01s 

ars  =  [36000, 36000,  -300,  -300] 
xsh  =  [36100, 36100,  -310,  -310] 
phi  =  0.2rad/s 


sample  time 
measurement  time 
initial  target  state 
initial  target  state  estimate 
turn  rate 


Aa  — 


mode  transition  probability  matrix  (3x3). 


If  in  straight  flight  mode: 


0.98 

0.01 


straight  path  will  continue 
change  to  left  (or  right)  turn  mode 


If  in  a  turn  flight  mode: 


0.99 

0.01 

Al ,  As  and  Ar= 

As 

Al 

At 

fh 


that  the  turn  will  continue 
change  to  straight  flight  mode. 

Probability  matrix  for  manoeuvre  induced  bin 

number  change,  each  (12  x  12). 

straight  path  gives  no  change  (Identity  matrix). 

left  turn  -  probability  of  transition  to  next  left  bin  is  0.05. 

right  turn  -  probability  of  transition  to  next  right  bin  is  0.05 

bin  transition  indicator  (36x1). 

Initialised  so  that  fh( 8)  =  1,  fh(i)  —  0  for  all  other  i 


E.2  Filters  in  Guidance  Loop 

h  =  0.0001 
hm  —  h 
radarh  =  1 
controlh  —  0.01 
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phi  =  0.2  rad/s 

war  =  (40m)2  target  range  measurement  variance 

ivar  =  (1.75  *  10  3ra<f/s)2  target  bearing  measurement  variance 

timesieps  =  100, 000 

ss  =  [2000, 1000, 660  x  cos(225),660  x  sm(225)J 

[100, 100, 10, 10]  initial  state  estimate  errors 

ia  —  225°  initial  target  heading 

xis  =  [0,0, 1000  x  cos(75),  1000  x  sin{ 75)] 

tiO  =  75°  interceptor  initial  bearing 
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valued  states)  and  considers  new  associated  approaches  such  as  the  polymorphic  estimator. 

Although  simulation  studies  were  performed,  the  polymorphic  estimator  had  serious  numeric  problems 
that  suggested  the  estimator  should  not  be  used  until  the  approach  is  refined.  This  report  is  intended  to 
facilitate  further  discussion  and  development  of  the  approach  (if  deemed  necessary);  hence,  it  includes 

details  of  the  assumptions  made  and  the  Matlab™  implementation. 

Page  classification:  UNCLASSIFIED 


